#pragma once

#include "pbrt.h"

namespace Pbrt
{
    template<typename T>
    inline bool isNaN(const T x)
    {
        return std::isnan(x);
        // For IEEE 754 float number : 
        // Check is NaN : x != x
    }
    template<>
    inline bool isNaN(const int x)
    {
        return false;
    }

    // Vector Declarations
    template<typename T>
    class Vector2
    {
        static_assert(
            std::is_same_v<T, float> || std::is_same_v<T, double> || std::is_same_v<T, int>,
            "Vector2 only supports float and double and int types.");
        public:
           T x, y;
           Vector2() : x(0), y(0) {}
           Vector2(T xx, T yy) : x(xx), y(yy)
           {
               DCHECK(!HasNaNs());
           } 
           bool HasNaNs() const
           {
               return isNaN(x) || isNaN(y);
           }
           // To do -> Complete Point2 class
           // explicit Vector2(const Point2<T> &p);
           // explicit Vector2(const Point3<T> &p);
           // ------- Realease Type -> No Debug -------
           #ifndef NDEBUG
           Vector2(const Vector2<T>& v)
           {
               DCHECK(!v.HasNaNs());
               x = v.x;
               y = v.y;
           }
           Vector2<T>& operator=(const Vector2<T>& v)
           {
               DCHECK(!v.HasNaNs());
               x = v.x;
               y = v.y;
               return *this;
           }
           #endif
           // --------------

           Vector2<T> operator+=(const Vector2<T>& v)
           {
                DCHECK(!v.HasNaNs());
                x += v.x;
                y += v.y;
                return *this;
           }
           Vector2<T> operator+(const Vector2<T>& v) const
           {
                DCHECK(!v.HasNaNs());
                return Vector2<T>(x + v.x, y + v.y);
           }

           Vector2<T> operator-(const Vector2<T>& v) const
           {
                DCHECK(!v.HasNaNs());
                return Vector2<T>(x - v.x, y - v.y);
           }
           Vector2<T> operator-=(const Vector2<T>& v)
           {
                DCHECK(!v.HasNaNs());
                x -= v.x;
                y -= v.y;
                return *this;
           }

           bool operator==(const Vector2<T>& v) const
           {
               return x == v.x && y == v.y;
           }
           bool operator!=(const Vector2<T>& v) const
           {
               return x != v.x || y != v.y;
           }

            template<typename U>
            Vector2<T> operator*(U f) const
            {
                DCHECK(!isNaN(f));
                return Vector2<T>(x * f, y * f);
            }

            template<typename U>
            Vector2<T>& operator*=(U f)
            {
                DCHECK(!isNaN(f));
                x *= f;
                y *= f;
                return *this;
            }

            template<typename U>
            Vector2<T> operator/(U f) const
            {
                CHECK_NE(f, 0);
                float inv = static_cast<float>(1.0 / f);
                return Vector2<T>(x * inv, y * inv);
            }

            template<typename U>
            Vector2<T>& operator/=(U f)
            {
                CHECK_NE(f, 0);
                float inv = static_cast<float>(1.0f);
                x *= inv;
                y *= inv;
                return *this;
            }

            Vector2<T> operator-() const
            {
                return Vector2<T>(-x, -y);
            }

            T operator[](std::size_t i) const
            {
                DCHECK(i >= 0 && i <= 1);
                if(i == 0) return x;
                return y;
            }

            T& operator[](std::size_t i)
            {
                DCHECK(i >= 0 && i <= 1);
                if(i == 0) return x;
                return y;
            }

            template<typename U>
            explicit operator Vector2<U>() const
            {
                return Vector2<U>(x, y);
            }

            float LengthSquared() const
            {
                return x * x + y * y;
            }

            float Length() const 
            {
                return std::sqrt(LengthSquared());
            }
    };

    template<typename T>
    inline std::ostream& operator<<(std::ostream& os, const Vector2<T>& vec)
    {
        os << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Vector2<float>& vec)
    {
        os << std::fixed << std::setprecision(9) << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Vector2<double>& vec)
    {
        os << std::fixed << std::setprecision(17) << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<typename T>
    class Vector3
    {
        static_assert(
            std::is_same_v<T, float> || std::is_same_v<T, double> || std::is_same_v<T, int> ,
            "Vector3 only supports float, double and int"
        );

        public:
            T x, y, z;

            Vector3() : x(0), y(0), z(0) {}
            Vector3(T x, T y, T z) : x(x), y(y), z(z) 
            {
                DCHECK(!HasNaNs());
            }
            bool HasNaNs() const
            {
                return isNaN(x) || isNaN(y) || isNaN(z);
            }
            // To do!
            // explicit Vector3(const Point3<T>& p);

            // Not define NoDebug -> Debug Type 
#ifndef NDEBUG
            Vector3(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
            }

            Vector3<T>& operator=(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
                return *this;
            }
#endif
            Vector3<T> operator+(const Vector3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Vector3<T>(x + v.x, y + v.y, z + v.z);
            }

            Vector3<T>& operator+=(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x += v.x;
                y += v.y;
                z += v.z;
                return *this;
            }

            Vector3<T> operator-(const Vector3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Vector3<T>(x - v.x, y - v.y, z - v.z);
            }
            
            Vector3<T>& operator-=(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x -= v.x;
                y -= v.y;
                z -= v.z;
                return *this;
            }

            bool operator==(const Vector3<T>& v) const
            {
                return x == v.x && y == v.y && z == v.z;
            }

            bool operator!=(const Vector3<T>& v) const
            {
                return !(*this == v);
            }

            template<typename U>
            Vector3<T> operator*(U s) const
            {
                return Vector3<T>(s * x, s * y, s * z);
            }

            template<typename U>
            Vector3<T>& operator*=(U s)
            {
                DCHECK(!isNaN(s));
                x *= s;
                y *= s;
                z *= s;
                return *this;
            }

            template<typename U>
            Vector3<T> operator/(U f) const
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                return Vector3<T>(inv * x, inv * y, inv * z);
            }

            template<typename U>
            Vector3<T>& operator/=(U f)
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                x *= inv;
                y *= inv;
                z *= inv;
                return *this;
            }

            Vector3<T> operator-() const
            {
                return Vector3<T>(-x, -y, -z);
            }

            float LengthSquared() const
            {
                return x * x + y * y + z * z;
            }

            float Length() const
            {
                return std::sqrt(LengthSquared());
            }

            // To do!
            // explicit Vector3(const Normal<T>& n);
    };

    template<typename T>
    inline std::ostream& operator<<(std::ostream& os, const Vector3<T>& v)
    {
        os << "[ " << v.x << ", " << v.y << ", " << v.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Vector3<float>& vec)
    {
        os << std::fixed << std::setprecision(9);
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Vector3<double>& vec)
    {
        os << std::fixed << std::setprecision(17); 
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z <<" ]";
        return os;
    }

    using Vector2f = Vector2<float>;
    using Vector3f = Vector3<float>;
    using Vector2i = Vector2<int>;
    using Vector3i = Vector3<int>;

    template<typename T>
    class Point2
    {
        public:
            T x, y;

            // Check Method
            bool HasNaNs() const
            {
                return isNaN(x) || isNaN(y);
            }

            // Constructor
            explicit Point2(const Point3<T>& p) : x(p.x), y(p.y)
            {
                DCHECK(!HasNaNs());
            }
            Point2() : x(0), y(0) {}
            Point2(T xx, T yy) : x(xx), y(yy)
            {
                DCHECK(!HasNaNs());
            }

            // constructor from point2<double> to point2<float>
            template<typename U>
            explicit Point2(const Point2<U>& p)
            {
                x = static_cast<T>(p.x);
                y = static_cast<T>(p.y);
                DCHECK(!HasNaNs());
            }

            // constructor from vector2<double> to vector2<float>
            template<typename U>
            explicit Point2(const Vector2<U>& p)
            {
                x = static_cast<T>(p.x);
                y = static_cast<T>(p.y);
                DCHECK(!HasNaNs());
            }

            // static_cast from point2<double> to vector2<double>
            template<typename U>
            explicit operator Vector2<U>() const
            {
                return Vector2<U>(
                    static_cast<U>(x),
                    static_cast<U>(y)
                );
            }

#ifndef NDEBUG
            Point2(const Point2<T>& p)
            {
                DCHECK(!p.HasNaNs());
                x = p.x;
                y = p.y;
            }

            Point2<T>& operator=(const Point2<T>& p)
            {
                DCHECK(!p.HasNaNs());
                x = p.x;
                y = p.y;
                return *this;
            }
#endif
            
            // point2 + vector2 : ray : o + td
            // o is orgin
            // t is time
            // d is unit direction
            // o + td is the ray destination
            Point2<T> operator+(const Vector2<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Point2<T>(x + v.x, y + v.y);
            }

            Point2<T>& operator+=(const Vector2<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x += v.x;
                y += v.y;
                return *this;
            }

            Point2<T> operator+(const Point2<T>& p) const
            {
                DCHECK(!p.HasNaNs());
                return Point2<T>(x + p.x, y + p.y);
            }
            Point2<T>& operator+=(const Point2<T>& p)
            {
                DCHECK(!p.HasNaNs());
                x += p.x;
                y += p.y;
                return *this;
            }

            // point2 - point2 -> represent a vector
            Vector2<T> operator-(const Point2<T>& p) const
            {
                DCHECK(!p.HasNaNs());
                return Vector2<T>(x - p.x, y -p.y);
            }

            // point2 - vector2 -> same as point2 + (-vector2)
            Point2<T> operator-(const Vector2<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Point2<T>(x - v.x, y - v.y);
            }

            Point2<T>& operator-=(const Vector2<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x -= v.x;
                y -= v.y;
                return *this;
            }

            Point2<T> operator-() const
            {
                return Point2<T>(-x, -y);
            }

            template<typename U>
            Point2<T> operator*(U f) const
            {
                return Point2<T>(f * x, f * y);
            }

            template<typename U>
            Point2<T>& operator*=(U f)
            {
                x *= f;
                y *= f;
                return *this;
            }

            template<typename U>
            Point2<T> operator/(U f) const
            {
                CHECK_NE(f, 0);
                float inv = static_cast<float>(1.0 / f);
                return Point2<T>(x * inv, y * inv);
            }

            template<typename U>
            Point2<T>& operator/=(U f)
            {
                CHECK_NE(f, 0);
                float inv = static_cast<float>(1.0 / f);
                x *= inv;
                y *= inv;
                return *this;
            }

            T operator[](std::size_t i) const
            {
                DCHECK(i >= 0 && i <= 1);
                if(i == 0) return x;
                return y;
            }

            T& operator[](std::size_t i)
            {
                DCHECK(i >= 0 && i <= 1);
                if(i == 0) return x;
                return y;
            }

            bool operator==(const Point2<T>& p) const
            {
                return x == p.x && y == p.y;
            }

            bool operator!=(const Point2<T>& p) const
            {
                return x != p.x || y != p.y;
            }
    };

    template<typename T>
    inline std::ostream& operator<<(std::ostream& os, const Point2<T>& vec)
    {
        os << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Point2<float>& vec)
    {
        os << std::fixed << std::setprecision(9) << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Point2<double>& vec)
    {
        os << std::fixed << std::setprecision(17) << "[ " << vec.x << ", " << vec.y << " ]";
        return os;
    }

    template<typename T>
    class Point3
    {
        static_assert(
            std::is_same_v<T, float> || std::is_same_v<T, double> || std::is_same_v<T, int> ,
            "Vector3 only supports float, double and int"
        );

        public:
            T x, y, z;

            Point3() : x(0), y(0), z(0) {}
            Point3(T x, T y, T z) : x(x), y(y), z(z) 
            {
                DCHECK(!HasNaNs());
            }
            bool HasNaNs() const
            {
                return isNaN(x) || isNaN(y) || isNaN(z);
            }
            // Point3<double> construct Point3<float>
            template<typename U>
            explicit Point3(
                const Point3<U>& p
            ) : x(static_cast<T>(p.x)),
                y(static_cast<T>(p.y)),
                z(static_cast<T>(p.z))
            {
                DCHECK(!HasNaNs());
            }
            // Convert from point3 to vector3
            template<typename U>
            explicit operator Vector3<U>() const
            {
                return Vector3<U>(x, y, z);
            }
#ifndef NDEBUG
            Point3(const Point3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
            }

            Point3<T>& operator=(const Point3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
                return *this;
            }
#endif
            // point3 + vector3 -> point3
            Point3<T> operator+(const Vector3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Point3<T>(x + v.x, y + v.y, z + v.z);
            }

            // point3 + vector3 -> point3
            Point3<T>& operator+=(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x += v.x;
                y += v.y;
                z += v.z;
                return *this;
            }

            // point3 - point3 -> vector3
            Vector3<T> operator-(const Point3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Point3<T>(x - v.x, y - v.y, z - v.z);
            }
            
            // point3 - vector3 -> new point3
            Point3<T> operator-(const Vector3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Point3<T>(x - v.x, y - v.y, z - v.z);
            }

            // point3 - vector3 -> this point3
            Point3<T>& operator-=(const Vector3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x -= v.x;
                y -= v.y;
                z -= v.z;
                return *this;
            }

            // point3 + point3 -> new point3
            Point3<T> operator+(const Point3<T>& p) const
            {
                DCHECK(!p.HasNaNs());
                return Point3<T>(x + p.x, y + p.y, z + p.z);
            }

            // point3 + point3 -> this point3
            Point3<T>& operator+=(const Point3<T>& p)
            {
                DCHECK(!p.HasNaNs());
                x += p.x;
                y += p.y;
                z += p.z;
                return *this;
            }

            template<typename U>
            Point3<T> operator*(U s) const
            {
                return Point3<T>(s * x, s * y, s * z);
            }

            template<typename U>
            Point3<T>& operator*=(U s)
            {
                DCHECK(!isNaN(s));
                x *= s;
                y *= s;
                z *= s;
                return *this;
            }

            template<typename U>
            Point3<T> operator/(U f) const
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                return Point3<T>(inv * x, inv * y, inv * z);
            }

            template<typename U>
            Point3<T>& operator/=(U f)
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                x *= inv;
                y *= inv;
                z *= inv;
                return *this;
            }

            // For point No operator-() const

            T operator[](std::size_t i) const
            {
                DCHECK(i >= 0 && i <= 2);
                if(i == 0) return x;
                if(i == 1) return y;
                return z;
            }

            T& operator[](std::size_t i)
            {
                DCHECK(i >= 0 && i <= 2);
                if(i == 0) return x;
                if(i == 1) return y;
                return z;
            }

            bool operator==(const Point3<T>& v) const
            {
                return x == v.x && y == v.y && z == v.z;
            }

            bool operator!=(const Point3<T>& v) const
            {
                return !(*this == v);
            }
    };

    template<typename T>
    inline std::ostream& operator<<(std::ostream& os, const Point3<T>& v)
    {
        os << "[ " << v.x << ", " << v.y << ", " << v.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Point3<float>& vec)
    {
        os << std::fixed << std::setprecision(9);
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Point3<double>& vec)
    {
        os << std::fixed << std::setprecision(17); 
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z <<" ]";
        return os;
    }

    using Point2f = Point2<float>;
    using Point2i = Point2<int>;
    using Point3f = Point3<float>;
    using Point3i = Point3<int>;

    template<typename T>
    class Normal3
    {
        static_assert(
            std::is_same_v<T, float> || std::is_same_v<T, double> || std::is_same_v<T, int> ,
            "Normal3 only supports float, double and int"
        );

        public:
            T x, y, z;

            Normal3() : x(0), y(0), z(0) {}
            Normal3(T x, T y, T z) : x(x), y(y), z(z) 
            {
                DCHECK(!HasNaNs());
            }
            bool HasNaNs() const
            {
                return isNaN(x) || isNaN(y) || isNaN(z);
            }

#ifndef NDEBUG
            Normal3(const Normal3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
            }

            Normal3<T>& operator=(const Normal3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x = v.x;
                y = v.y;
                z = v.z;
                return *this;
            }
#endif
            // construct from vector3
            explicit Normal3<T>(const Vector3<T>& v)
            : x(v.x), y(v.y), z(v.z)
            {
                DCHECK(!v.HasNaNs());
            }

            // normal3 + normal3 -> new normal3
            Normal3<T> operator+(const Normal3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Normal3<T>(x + v.x, y + v.y, z + v.z);
            }

            // normal3 + normal3 -> this normal3
            Normal3<T>& operator+=(const Normal3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x += v.x;
                y += v.y;
                z += v.z;
                return *this;
            }

            // normal3 - normal3 -> new normal3
            Normal3<T> operator-(const Normal3<T>& v) const
            {
                DCHECK(!v.HasNaNs());
                return Normal3<T>(x - v.x, y - v.y, z - v.z);
            }
            
            // normal3 - normal3 -> this normal3
            Normal3<T>& operator-=(const Normal3<T>& v)
            {
                DCHECK(!v.HasNaNs());
                x -= v.x;
                y -= v.y;
                z -= v.z;
                return *this;
            }

            bool operator==(const Normal3<T>& v) const
            {
                return x == v.x && y == v.y && z == v.z;
            }

            bool operator!=(const Normal3<T>& v) const
            {
                return !(*this == v);
            }

            // scale normal -> new normal
            template<typename U>
            Normal3<T> operator*(U s) const
            {
                return Normal3<T>(s * x, s * y, s * z);
            }

            // scale normal -> this normal
            template<typename U>
            Normal3<T>& operator*=(U s)
            {
                DCHECK(!isNaN(s));
                x *= s;
                y *= s;
                z *= s;
                return *this;
            }

            // normal / scale -> new normal
            template<typename U>
            Normal3<T> operator/(U f) const
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                return Normal3<T>(inv * x, inv * y, inv * z);
            }

            // normal / scale -> this normal
            template<typename U>
            Normal3<T>& operator/=(U f)
            {
                CHECK_NE(f, 0);
                float inv = 1.0f / f;
                x *= inv;
                y *= inv;
                z *= inv;
                return *this;
            }

            float LengthSquared() const
            {
                return x * x + y * y + z * z;
            }

            float Length() const
            {
                return std::sqrt(LengthSquared());
            }

            T operator[](std::size_t i) const
            {
                DCHECK(i >= 0 && i < 2);
                if(i == 0) return x;
                if(i == 1) return y;
                else return z;
            }

            T& operator[](std::size_t i)
            {
                DCHECK(i >= 0 && i < 2);
                if(i == 0) return x;
                if(i == 1) return y;
                else return z;
            }

    };

    template<typename T>
    inline std::ostream& operator<<(std::ostream& os, const Normal3<T>& v)
    {
        os << "[ " << v.x << ", " << v.y << ", " << v.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Normal3<float>& vec)
    {
        os << std::fixed << std::setprecision(9);
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z << " ]";
        return os;
    }

    template<>
    inline std::ostream& operator<<(std::ostream& os, const Normal3<double>& vec)
    {
        os << std::fixed << std::setprecision(17); 
        os << "[ " << vec.x << ", " << vec.y << ", " << vec.z <<" ]";
        return os;
    }

    using Normal3f = Normal3<float>;
    // -----------------------------------

    // Axis Aligned Bounding Box (Two Dimensions)
    template<typename T>
    class Bounds2
    {
        public:
            Point2<T> min, max;
            Bounds2()
            {
                T minNum = std::numeric_limits<T>::min();
                T maxNum = std::numeric_limits<T>::max();
                min = Point2<T>(maxNum, maxNum);
                max = Point2<T>(minNum, minNum);
            }

            // A point is a minimal bounding box 2 dimensions
            explicit Bounds2(const Point2<T>& p)
            : min(p), max(p) {}

            Bounds2(const Point2<T>& p1, const Point2<T>& p2)
            {
                min = Point2<T>(std::min(p1.x, p2.x), std::min(p1.y, p2.y));
                max = Point2<T>(std::max(p1.x, p2.x), std::max(p1.y, p2.y));
            }

            // tranform Bounds2<float> to Bounds2<double>
            template<typename U>
            explicit operator Bounds2<U>() const
            {
                return Bounds2<U>(
                    static_cast<Point2<U>>(min),
                    static_cast<Point2<U>>(max)
                );
            }

            // return a vector point from min point to max point
            Vector2<T> Diagonal() const
            {
                return max - min;
            }

            T Area() const
            {
                auto dia = Diagonal();
                return dia.x * dia.y;
            }

            int Max_Extend_Axis() const
            {
                auto diagonal = Diagonal();
                if(diagonal.x > diagonal.y) return 0;
                else return 1;
            }

            // For none normal type
            // const call operator[]
            // always return a const reference
            inline const Point2<T>& operator[](
                int i
            ) const
            {
                DCHECK(i == 0 || i == 1);
                return i == 0 ? min : max;
            }

            // None const call operator[]
            // return a reference
            inline Point2<T>& operator[](
                int i
            )
            {
                DCHECK(i == 0 || i == 1);
                return i == 0 ? min : max;
            }

            bool operator==(const Bounds2<T>& b) const
            {
                return b.min == min && b.max == max;
            }

            bool operator!=(const Bounds2<T>& b) const
            {
                return !(*this == b);
            }

            Point2<T> Lerp(const Point2f& t) const
            {
                return Point2<T>(
                    Lerp(t.x, min.x, max.x),
                    Lerp(t.y, min.y, max.y)
                );
            }

            Point2<T> Offset(const Point2<T>& p) const
            {
                Vector2<T> o = p - min;
                if(max.x > min.x) o.x /= max.x - min.x;
                if(max.y > min.y) o.y /= max.y - min.y;
                return o;
            }

            void BoundSphere(Point2<T>& center, float& rad) const
            {
                center = (min + max) / 2;
                rad = Inside(center, *this) ? Distance(center, min) : 0;
            }

            friend std::ostream& operator<<(std::ostream& os, const Bounds2<T>& b)
            {
                // Depend on type T to call different version of operator<< of point2
                os << "[" << b.min << "-" << b.max << "]";
                return os;
            }

    };

    template<typename T>
    class Bounds3
    {
        public:
            Point3<T> min, max;
            Bounds3()
            {
                T minNum = std::numeric_limits<T>::min();
                T maxNum = std::numeric_limits<T>::max();
                min = Point3<T>(maxNum, maxNum, maxNum);
                max = Point3<T>(minNum, minNum, minNum);
            }

            // A point is a minimal bounding box 3 dimensions
            explicit Bounds3(const Point3<T>& p)
            : min(p), max(p) {}

            // Two points define a Bounding Box
            Bounds3(const Point3<T>& p1, const Point3<T>& p2)
            {
                min = Point3<T>(
                    std::min(p1.x, p2.x), std::min(p1.y, p2.y), std::min(p1.z, p2.z)
                );
                max = Point3<T>(
                    std::max(p1.x, p2.x), std::max(p1.y, p2.y), std::max(p1.z, p2.z)
                );
            }

            // tranform Bounds3<float> to Bounds3<double>
            template<typename U>
            explicit operator Bounds3<U>() const
            {
                return Bounds3<U>(
                    static_cast<Point3<U>>(min),
                    static_cast<Point3<U>>(max)
                );
            }

            // Get Bounds3 Corner coordinates
            Point3<T> Corner(int corner) const
            {
                DCHECK(corner >= 0 && corner <= 7);
                return Point3<T>((*this)[(corner & 1)].x,
                         (*this)[(corner & 2) ? 1 : 0].y,
                         (*this)[(corner & 4) ? 1 : 0].z);
            }

            // Diagonal Vector
            Vector3<T> Diagonal() const
            {
                return max - min;
            }

            T Area() const
            {
                Vector3<T> d = Diagonal();
                return 2 * (d.x * d.y + d.x * d.z + d.y * d.z);
            }

            T Volume() const
            {
                Vector3<T> d = Diagonal();
                return d.x * d.y * d.z;
            }

            int Max_Extend_Axis() const
            {
                Vector3<T> diagonal = Diagonal();
                if(diagonal.x > diagonal.y && diagonal.x > diagonal.z)
                    return 0;
                else if(diagonal.y > diagonal.z)
                    return 1;
                else
                    return 2;
            }

            inline const Point3<T>& operator[](int i) const
            {
                DCHECK(i == 0 || i == 1);
                return i == 0 ? min : max;
            }

            inline Point3<T>& operator[](int i)
            {
                DCHECK(i == 0 || i == 1);
                return i == 0 ? min : max;
            }

            bool operator==(const Bounds3<T>& b) const
            {
                return min == b.min && max == b.max;
            }

            bool operator!=(const Bounds3<T>& b) const
            {
                return min != b.min || max != b.max;
            }

            Point3<T> Lerp(const Point3f& t) const
            {
                return Point3<T>(
                    Lerp(t.x, min.x, max.x),
                    Lerp(t.y, min.y, max.y),
                    Lerp(t.z, min.z, max.z)
                );
            }

            Point3<T> Offset(const Point3<T>& p) const
            {
                Vector3<T> o = p - min;
                if(max.x > min.x)
                    o.x /= max.x - min.x;
                if(max.y > min.y)
                    o.y /= max.y - min.y;
                if(max.z > min.z)
                    o.z /= max.z - min.z;
                return o;
            }

            void BoundSphere(Point3<T>& center, float& rad) const
            {
                center = (min + max) / 2;
                rad = Inside(center, *this) ? Distance(center, min) : 0;
            }

            friend std::ostream& operator<<(std::ostream& os, const Bounds3<T>& b)
            {
                os << "[" << b.min << "-" << b.max << "]";
                return os;
            }
            
            // Wait to Implementation
            // when the Ray class is implementation
            inline bool IntersectP(
                const Ray& ray, float& hit0, float hit1&
            ) const;

            inline bool IntersectP(
                const Ray& ray, const Vector3f& invDir, const int dirIsNeg[3]
            ) const;
    };

    using Bounds2f = Bounds2<float>;
    using Bounds3f = Bounds3<float>;
    using Bounds2i = Bounds2<int>;
    using Bounds3i = Bounds3<int>;

    // Just a tag
    // Not follow the rule of stl standard
    class Bounds2iIterator : public std::forward_iterator_tag
    {
        private:
            Point2i p;
            const Bounds2i* bounds;

            // Iterator over the Int bound 
            void advance()
            {
                ++p.x;
                if(p.x == bounds->max.x)
                {
                    p.x = bounds->min.x;
                    ++p.y;
                }
            }

        public:
            Bounds2iIterator(
                const Bounds2i& b,
                const Point2i& p
            ) : p(p), bounds(&b) {}

            // pre increment
            Bounds2iIterator& operator++()
            {
                advance();
                return *this;
            }

            Bounds2iIterator operator++(int)
            {
                auto temp = *this;
                advance();
                return temp;
            }

            bool operator==(const Bounds2iIterator& it) const
            {
                return p == it.p && bounds == it.bounds;
            }

            bool operator!=(const Bounds2iIterator& it) const
            {
                return p != it.p || bounds != it.bounds;
            }

            Point2i operator*() const
            {
                return p;
            }

    };

    template<typename T>
    bool Inside(const Point2<T>& p, const Bounds2<T>& b)
    {
        return p.x >= b.min.x && p.x <= b.max.x &&
               p.y >= b.min.y && p.y <= b.max.y;
    }

    template<typename T>
    bool Inside(const Point3<T>& p, const Bounds3<T>& b)
    {
        return p.x >= b.min.x && p.x <= b.max.x &&
               p.y >= b.min.y && p.y <= b.max.y &&
               p.z >= b.min.z && p.z <= b.max.z;
    }

    template<typename T>
    float Distance(const Point2<T>& p1, const Point2<T>& p2)
    {
        return (p1 - p2).length();
    }

    template<typename T>
    float Distance(const Point3<T>& p1, const Point3<T>& p2)
    {
        return (p1 - p2).length();
    }



};